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Abstract 

The hydrodynamic equations for a model of a confined quasi-two-dimensional gas of smooth in¬ 
elastic hard spheres are derived from the Boltzmann equation for the model, using a generalization 
of the Chapman-Enskog method. The heat and momentum fluxes are calculated to Navier-Stokes 
order, and the associated transport coefficients are explicitly determined as functions of the coeffi¬ 
cient of normal restitution and the velocity parameter involved in the definition of the model. Also 
an Euler transport term contributing to the energy transport equation is considered. This term 
arises from the gradient expansion of the rate of change of the temperature due to the inelasticity 
of collisions, and vanishes for elastic systems. The hydrodynamic equations are particularized for 
the relevant case of a system in the homogeneous steady state. The relationship with previous 
works is analyzed. 

PACS numbers: 45.70.Mg,45.70.-n,05.20.Jj,51.10,+y 
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I. INTRODUCTION 


Granular gases frequently exhibit flows similar to those of normal gases, and for practical 


purposes t 
equations 


res^ flows are often successfully described by phenomenological hydrodynamic 
1 


, 2], The basis for such macroscopic balance equations are in the more fun¬ 
damental statistical mechanics and kinetic theory descriptions of granular gases. In this 
context, the idealized model of a granular gas as a monodisperse system of smooth inelastic 
hard spheres or disks with a constant coefficient of normal restitution has been extensively 


employed 


, 4]. For this system, hydrodynamic equations to Navier-Stokes order have been 


derived with expressions for the parameters appearing in them. Starting from the Boltzmann 


equation for inelastic hard spheres or disks 


, 6] and also from the revised Enskog theory 


7), the transport coefficients have been evaluated by using an extension of the Chapman- 
Enskog method. The predictions from the Boltzmann equation have been found to be in 
good agreement with the values obtained by particle simulation methods in the dilute limit 
6]. Using linear response theory, formal Green-Kubo like expressions for the transport 
coefficient^have been derived for low density granular gases 8, 9J, and also for arbitrary 
densities [lo|, ll|. The latter are not tied to any specific kinetic equation, but their explicit 
evaluation requires the introduction of some approximations [loj |. 

In normal fluids, non-equilibrium steady states can be generated by imposing appropriate 
boundary conditions. Moreover, the control of the boundary conditions permits to keep the 
gradients of the hydrodynamic fields small, so that the steady state can be studied in the 
Navier-Stokes domain of the hydrodynamic equations. In granular gases, a new class of 
steady states shows up. In them, stationarity is reached by an autonomous balance between 
external constrains and the internal cooling. A typical example is a system under shear 
flow. There is viscous heating due to the work made on the system at the boundaries. If the 
system is a granular gas, stationarity is possible when the viscous heating is compensated 
by the dissipation due to the inelasticity of collisions. The steady state has uniform density 


and temperature, and a flow ve 
it has been extensively studied 


ocity with a linear profile. Due to its macroscopic simplicity, 


12-[l6j. This particular steady state exemplifies two features 


that are characteristic of hydrodynamic steady states of granular gases. First, to compensate 
the energy dissipation in collisions, the system must develop spatial gradients generating an 
energy flux, i.e. it must be inhomogeneous. Second, the energy balance leads to a coupling 
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between gradients and inelasticity, so that the limit of small gradients also implies the quasi¬ 
elastic limit. Even more, the above coupling is often non-analytical [blj], implying that the 
macroscopic description of the steady state can never be brought within the range of validity 
of the Navier-Stokes hydrodynamics in those cases. 

An interesting alternative to the kind of steady statesjaf granular gases described above 


has being attracting increasing interest in the last years 


3 - 21 ]. 


A granular gas is confined 


to a quasi-two-dimensional geometry by placing it between two large parallel plates in the 
horizontal directions, while the distance between the two plates is smaller than two particle 
diameters, so the system is actually a monolayer since the particles can not jump on another. 
The container is vertically vibrated to inject energy through the collisions of the particles 
with the top and bottom walls. The two-dimensional dynamics of the system when seen from 
above is considered. It has been observed that it corresponds to that of a two-dimensional 
granular fluid. Moreover, the system remains homogeneous under a large range of parameters 
and it eventually reaches a steady homogeneous state. Very recently 22], an idealized model 
has been proposed trying to describe the horizontal dynamics in the above experiment, 
assuming that the particles are smooth inelastic hard spheres. Then, the projections of the 
particles on the horizontal plane are described as inelastic hard disks, whose collision rule 
is modified in order to incorporate a mechanism to transfer the energy injected vertically to 
the horizontal degrees of freedom. In this sense, it can be classified as a collisional model. 
The new collision rule has a constant coefficient of normal restitution a, and contains a 
characteristic velocity A that is added to each particle in a collision so that the normal 
component of the relative velocity is increased by 2A in the collision, independently of 
the effect of the coefficient of normal restitution. The methods of non-equilibrium statistical 


mechanics and kinetic t 


applied to the model 


reory developed for inelastic hard spheres and disks 


23 


24] have been 


. In this way, the Boltzmann equation and the revised Enskog theory 


have been formulated. Moreover, the existence of homogeneous hydrodynamics has been 
analyzed 26]. There is a time regime over which the granular temperature of a homogeneous 
system obeys a closed (hydrodynamic) equation. In the long time limit, the temperature 
tends to its steady value. Moreover, it has been shown that the homogeneous relaxation 
of the temperature of the system presents nonlinear memory effects 27], which can be 
considered as reminiscent of the Kovacs memory effect occurring in the relaxation towards 
equilibrium of molecular fluids 28]. 
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The aim of this paper is to derive the hydrodynamic equations to Navier-Stokes order for 
a dilute confined granular gas as described by the collisional model proposed by Brito et al 
22]. The starting point will be the Boltzmann kinetic equation and the method to be used 
a generalization of the Chapman-Enskog procedure. The derivation is based on a special 
“normal” solution of the kinetic equation expanded to low order in the gradients of the 
hydrodynamic fields. The zeroth order approximation is not a local version, both in space 
and time, of the distribution function of the homogeneous steady state, but it is based on the 
distribution describing the homogeneous hydrodynamics. This is an important conceptual 
and practical difference with the standard application of the Chapman-Enskog method to 
molecular systems. Of course, it is also possible to consider states close to a stationary one 
and carry out, for instance, linear response theory around that state to compute transport 
properties associated to that particular state. The ranges of applicability of the results 
obtained by both methods are clearly different, although there can be a common limit for 
the simultaneous validity of both. In particular, the Navier-Stokes shear viscosity of a 
confined dilute granular gas described by the collisional model in a stationary and uniform 


Couette flow has been computed by employing linear response theory 29]. 
obtained here for the shear viscosity will be related with those reported in Ref. 


he results 


29]. 


The remaining of this paper is organized as follows. In the next section, the Boltzmann 
equation for the model is given, and the exact balance equations for mass, momentum, 
and energy are derived from it. The Chapman-Enskog method for obtaining a “normal” 
solution of the kinetic equation as an expansion in the gradients of the hydrodynamic fields 
is described. Results through Navier-Stokes order for the pressure tensor and the heat flux 
are given. The associated transport coefficients are shown to obey a complete set of first 
order differential equations. Some details of the calculations are given in Appendices [A] and 
iBl while the explicit results for the transport coefficients are presented in Sec. IIIII The 
theory is not restricted to any range of values of the coefficient of normal restitution nor 
of the characteristic velocity of the model A. The contributions to the transport equations 
coming from the energy sink term due to the non conservation of kinetic energy in collisions 
are also discussed. It is shown that there is a first order in the gradients contribution, Euler 
term, which is proportional to the divergence of the velocity field. The associated transport 
coefficient is explicitly evaluated. Appendix [C] provides an sketch of the calculation of this 
coefficient. The Euler term does not exist in molecular systems and it is peculiar of inelastic 
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collisions 30 J 32|, although it vanishes in the low density limit of granular gases composed of 
smooth inelastic hard spheres or disks j^|. The expressions of the transport coefficients are 


particularized for the steady homogeneous state in Sec. IIII1 The peculiarity of 


which the temperature is a function of the parameters defining the system 


22 


his state, in 


25], leads to 


much simpler expressions of the coefficients. Further comments and conclusions, as well as 
the relationship with some previous work for the viscosity transport coefficient, are given in 
Sec. ED 

II. CHAPMAN ENSKOG SOLUTION OF THE BOLTZMANN EQUATION 

The Boltzmann equation obeyed by the one-particle distribution function, f(r,v,t ), of 
the model reads [25] 


d d 

dt + V ' dr )^( r, ' u,t ) = I dvi 


Here, Tq is the binary collision operator defined by 


( 1 ) 


Tq{v 1 , v 2 ) = a d 1 I da [0 (v l2 ■ a - 2 A) (v i2 ■ a - 2A) a 2 b a . 1 ( 1, 2) - 0(vi 2 • $){v 12 ■ S)\ , 

( 2 ) 

where cr is the diameter of the particles, d the dimension of the system 33J, a is the unit 


vector joining the center of the two particles at contact, V \ 2 = V\ — v 2 is the relative velocity, 
Q(x) is the Heaviside step function, a is the coefficient of normal restitution defined in the 
interval 0 < a < 1, A is some positive characteristic speed, and 6 (T 1 (f,j) is an operator 
changing all the velocities and Vj to its right into the precollisional values corresponding 
to a collision between them defined by er, i.e., 


*\ * 1 + a — . ^ 

b a ( i,j)Vi = v { = Vi --— Vij ■ a a +-, 

la a 


K (i,3)vj = = v j + 


1 ~t~ a 
2 a 


■ (T(T 


For arbitrary velocity functions, a(v iy Vj) and b(Vi,Vj), it is 25] 


Acf 


a 


( 3 ) 

( 4 ) 


dvi / dv j b(v i ,v j )T 0 (v i ,v j )a(v i ,v j ) = / dvi / dvj a(v h Vj)T 0 (vi, Vj)b(vi, Vj), (5) 


where 


T 0 (vi, vj) =a d 1 / da Q(-Vij ■ ct)| v {j ■ S\ j) - 1] 


( 6 ) 
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The operator b a (i,j) is the inverse of 6 <T 1 (f, j), i.e. it changes Vi and Vj into their postcolli- 
sional values, v[ and n', given by 


b<r (b j)Vi = v[ = Vi - 


1 + 01 - ^ 

—-—Vij ■ aa + A a 


b*(i,j)v j = = vj + 


1 T ex 


-v 




aa 


A a. 


The kinetic energy change in a collision is 

A 2 - a A v^ ■ a - (v^ - a) 2 , 

with m being the mass of a particle. Using the identity ([5]) it is easily found that 



(7) 

( 8 ) 

(9) 


dv / dv 1 T 0 (v,v 1 )f(r,v,t)f(r,v 1 ,t) = 0, 


( 10 ) 


dv / dv 1 vT 0 (v,v 1 )f(r,v,t)f(r,v 1 ,t) = 0, 


( 11 ) 


reflecting the conservation of the number of particles and the momentum, respectively. On 
the other hand, it is 


/ 771V _ 

dv / dv 1 —-T 0 (v,v 1 )f(r,v,t)f(r,v u t)=u{f,f]. 


( 12 ) 


The term uj[f, f] provides the rate of energy change due to the inelasticity of collisions, and 
the functional oj[f, h] is 


rS dj -\)/ 2 ma d ^ x f f 

uj[f,h\ = --- / dv i / dv 2 f(r,vi,t)h(r,v 2 ,t) 


x 


A 2 ui 2 7r 1 / 2 aAu 2 2 (1 — a 2 )vf 2 


(13) 


r (^±i) 2 r (^) 4 r (^) 

The macroscopic number of particles density, n(r,t), flow velocity, u(r,t), and granular 
temperature, T(r,t), are defined from the one-particle distribution function in the usual 
way, 

n(r,t) = [ dv f(r,v,t), (14) 


n(r,t)u(r,t) = J dvvf(r,v,t), 

d , , f , mV 2 . 

-n(r, t)T(r, t) = dv —- f(r, v, t), 


(15) 

(16) 
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where V(r,t) — v — u(r,t ) is the velocity of the particle relative to the flow field. Balance 
equations for the above fields follow by taking velocity moments in the Boltzmann equation, 
Eq. (P), 


dn 

~dt 


+ V • (riu) = 0, 


du 

~di 


u ■ V« + ( mn ) i V ■ P = 0, 


dT 

~di 


u • VT + — (P : Vu + V • Jo) = -C T. 
nd 


In the above equations, the pressure tensor, P, and the heat flux, J q , are defined by 


P (r,t)=m / dvV(r,t)V(r,t)f(r,v,t) 


(17) 

(18) 

(19) 

( 20 ) 


and 


777 / / 

J q{ r ^) = ~2 J ( 21 ) 

respectively. In addition, Eq. (1T9I) contains the rate of change of the temperature, £( r,t ), 
due to the inelasticity of collisions, whose expression is 

2 


C(r,t) = 


/]• 


( 22 ) 


n(r, t)T(r, t)d 

The minus sign has been introduced by analogy with a system of smooth inelastic hard 
spheres or disks, but in the present context it does not presuppose that ( is (semi)defined 
positive jl| - 

To close the balance equations (TTTTl-(TTUT) . it is necessary to express the fluxes and the 
temperature change rate in terms of the macroscopic fields, by means of some constitutive 


relations. To accomplish this, the Chapman-Enskog theory 


34] assumes the existence of a 


normal solution of the Boltzmann equation, i.e. a solution in which all the dependence of 
the distribution function on position and time occurs through its functional dependence on 
the macroscopic fields n, u , and T, 


f(r,v,t) = f\v\n,u,T}. 


(23) 


Next, it is assumed that the space and time variations of the fields are small, so that the 
functional dependence of the distribution function on the fields can be localized in space and 
time by means of an expansion in gradients. Then, the distribution function is expressed as 
a power series expansion in a formal uniformity parameter e, 

/ = / (0) + e/ {1) + e 2 / (2) + -- - . (24) 
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Since the aim is to generate a gradient expansion, a factor of e is assigned to every gradient 
operator. Moreover, the Chapman-Enskog method uses the multiple-scale perturbation 

into the definition 


theory 35]. In practice, this is done by using the expansion in Eq. 
of the fluxes and the dissipation rate (. Then the resulting expansions are introduced into 
the macroscopic balance equations to get an identification of the time derivatives of the 
macroscopic fields in the form of an expansion in the uniformity parameter, 


— = a (0) + e a (1) + e 2 a (2) + • • •. 


dt 


(25) 


Details of the application of the method are given in Appendices [A] [B] and [Q To first order 
in the gradients, the pressure tensor and heat flux are given by 


P = nT\ — r/ 


V-u + (Vu) + --V • u\ 


d 


(26) 

Jq = —kVT — pVn, (27) 

where I is the unit tensor in d dimensions, (Vu) + is the transposed of Vit, r] is the coefficient 
of shear viscosity, n the heat conductivity, and p a new coefficient coupling the heat flux and 


the density gradient, which is peculiar of inelastic collisions 


36]. To distinguish between the 


two energy transport coefficients, sometimes k is referred to as the thermal heat conductivity 
and /i as the diffusive heat conductivity. The transport coefficients are determined by the 
normal solutions of the first order differential equations 

d°)\ 2 5 / 2 7T^ 


C (0) A* di] 


dA* 


- _c 

Urj 


V = 


C (0) A* <97c 
2 dA* 


+ + 


A* <9C (0) 
~2 dA* 


-2C 


, 7 ( 0 ) 


K = 


(d + 2)T (d/2) ’ 

2 5 / 2 (d-l)vr i f i 
d(d+2)T(d/2) 


(28) 


l + 2<i2-4i|S: ], (29) 


2 dA* 


C (0) A* 

2 dA* 


+ V 


- 3C 


(o)' 


_ -(o)_ 2 5/2 (d - l) 7 r d 2 1 a 2 

H - C K = 


d(d + 2)T(d/2) ' 

In these equations, A* = A (m/2T) 1 ^ 2 and dimensionless transport coefficients have been 
introduced. They are defined by 


n = 


Vo 


K = 


where 


Vo = 


2 + d^ 


K 

«o’ 

■si 

III 

4 s 

I s ~ 

(31) 

d—1 

2 

(mT) 1/2 cr _ ( d_1 \ 

(32) 
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and 


did + 2) 2 

K o = 77m —7T r (d/2) 71 


d -1 
2 


rjl \ 1/2 

m ) 


a 


-id- 1 ) 


(33) 


16(d — 1) 

are the shear viscosity and the (thermal) heat conductivity, respectively, of a molecular gas 
described by the Boltzmann equation, with the Boltzmann constant set equal to unity. The 
dimensionless functions introduced in Eqs. (|28]) - (l30)) are 


c (0) = 


C c °) /mb / 2 _ 


Tiff 1 


d-i \2T7 


IS V = 


^7) 


(m ^V 2 _ _ _ v K (m^ 1 / 2 

u. = v„ = nad _ 1 ) 


(34) 


na d ~ 1 \2T 7 

The expression of the zeroth order rate of change of the temperature, is given in 

Eq. (IA12D . while the frequencies v v and u K are given in Eqs. (IB5D and ()B6j) . respectively. 
Some details of the calculations are shown in Appendices [A] and [0 Finally, Eqs. (l28lh 
(1301) have been obtained by considering the distribution function in the called first Sonine 
approximation, in which the deviation of the one-particle distribution function of the gas 
from the Gaussian are characterized by the coefficient ei 2 , given by the normal solution of 
the ordinary differential equation (1A13I) . Moreover, the same kind of approximation has 
been considering upon evaluating the hydrodynamic fluxes. This is the usual approximation 
to obtain explicit expressions for the transport coefficients of a gas with elastic collisions 
and there is no reason to question its accuracy here as well. Actually, it has been shown to 
lead to quantitatively right approximations in the case of a system of smooth inelastic hard 


spheres or disks 


6 ], 


III. EULER AND NAVIER-STOKES TRANSPORT COEFFICIENTS 

As a consequence of the confinement of the fluid and its description by means of a modified 
hard collision, there is a contribution to the hydrodynamic equations of the rate of change 
of the temperature £( r,t) of first order in the gradients, namely it is (see Eq. flA31jl ) 

C (1) M) = CiV-u. (35) 

The Euler transport coefficient £i represents dissipation due to the inelastic character of 
collisions proportional to V • u. It has no analogue for elastic fluids, where the Euler 
hydrodynamics (first order in the gradients of the fields) is referred to as the “perfect fluid” 
equations, since there is not dissipation in that limit. The expression derived here vanishes in 
the limit A* —» 0, as a consequence of symmetry considerations j^], i.e. there is no dissipation 


9 
















to Euler order in a dilute gas of smooth inelastic hard spheres or disks. On the other hand, a 
term of the form given in Eq. ([35]) is present in the hydrodynamic equations even in systems 
of smooth particles if density effects are considered [7, lll,l3Qj. The calculation of Ci using the 
Chapman-Enskog procedure, requires to determine The details of the calculation are 
given in Appendix O where the first Sonine polynomial expansion is again employed. The 
same approximation was used above to compute the Navier-Stokes transport coefficients. In 
Fig. [U the dimensionless coefficient Ci is plotted as a function of the speed parameter for 
two values of the coefficient of normal restitution, namely a = 0.8 and a = 0.9. The value 
of the transport coefficient in the homogenous steady state corresponding to each value of 
a is indicate in the figure. To fold the values of the transport coefficient a coupled of first 
order differential equations have to be solved, namely Eqs. (IA13D and (IC3j) . The curves 
reported in the figure correspond to the hydrodynamic solutions of the equations, i.e. they 
are identified independently of the initial conditions. The method is described below in 
detail for the shear viscosity coefficient and it will be not be discussed now. 

The coefficient is the contribution of the energy source in collisions to what would phys¬ 
ically constitute the effects of the hydrostatic pressure at the Euler order. If a small element 
of the confined granular gas is considered, then the pressure that the fluid element exerts on 
its boundaries is decreased or increased by the energy lost locally in collisions. At the level 
of linear hydrodynamics, the pressure and the Euler dissipative term are indistinguishable 
in their physical implications 11]. 

To compute the Navier-Stokes transport coefficients, Eqs. (j28l) - ([30l) have to be solved. In 
the equations, the hydrodynamic expression of the second Sonine coefficient, CZ 2 , has to be 
used. The latter is obtained by numerically solving Eq. (1A13[) as a function of A* for a fixed 
value of a, and a given initial condition 02 ( 0 , Aq) = 02 , 0 - It is seen that all the trajectories 
converge quite fast towards a universal curve, identified as the hydrodynamic expression of 
the second Sonine coefficient 26]. A similar method has been employed here to generate 
the hydrodynamic transport coefficients. Actually, what has been done is to simultaneously 
solve the equation for 02 and those for the transport coefficients. Since the rate of variation 
of the temperature vanishes in the steady state, the equations for the transport coefficients 
have a singularity at the steady value A* = A* t . Therefore, in the numerical simulations, 
trajectories have been generated starting from both Aq > A* t and with Ajj < A: t . The 


hydrodynamic solution giving the expression of the transport coefficient is the common 
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FIG. 1: (Color online) Dimensionless Euler transport coefficient (j as a function of dimensionless 
characteristic speed A* in a two-dimensional system. The (red) dashed line is for a coefficient of 
normal restitution a = 0.8 and the (black) solid line is for a = 0.9. The (blue) dots indicate the 
values of the transport coefficient in each of the two steady states. 

part of all the numerical solutions. As an example, the numerical results obtained for the 
dimensionless coefficient fj in a two dimensional system with a = 0.9 are shown in Fig. [21 
All the numerical trajectories converge towards the same curve, then forgetting the initial 
conditions used to generate them. This is consistent with the existence of a hydrodynamic 
shear viscosity being a function of only the local hydrodynamic fields, but not of the previous 
history or some initial values. In the particular case shown in Fig. [2j several initial values 
of the viscosity parameter corresponding to A* = 0.005 and A* = 10 have been employed. 
The curves tend to converge quite fast in the range 0 < A* < 0.2. For A* > 0.2, the 
dependence of the solution of the differential equation on the initial value of fj used for 
A* = Ag is rather strong and much more intensive numerical simulations would be needed 
to identify the value of the hydrodynamic shear viscosity. The two particular solutions 
drawn in Fig. [2] for A* > A* st correspond to fj(A* = 10) = 100 and fj(A* = 10) = 0, 
respectively, while the third plotted particular solution has been obtained with the initial 
condition rj( A* = 0.005) = 10. Results obtained with other initial conditions can not be 
distinguished from the normal curve on the scale of the figure. 

In Figs. [3] - [5] the coefficients of shear viscosity, fj, (thermal) heat conductivity, k, and 
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FIG. 2: (Color online) Adimensionalized quantity rj as a function of the dimensionless characteristic 
speed A* in a two-dimensional system with a = 0.9. The (red) dashed lines correspond to numerical 
solutions of Eq. (1281) obtained by using different initial conditions, i.e. different values for the pair 
Aq,? 7 (Aq). The (black) solid line is the universal curve to which all the solutions converge. This is 
precisely the function identified as the dimensionless hydrodynamic shear viscosity. Also indicated 
in the figure is the steady value of A*, denoted by A* t and the shear viscosity of the steady state, 
Vst ■ 


diffusive heat conductivity, /Z, are plotted as a function of the dimensionless characteristic 
speed for a two-dimensional system. Two values of the coefficient of normal restitution 
have been considered, namely a = 0.8 and a = 0.9. The reported curves correspond to 
the hydrodynamic transport coefficients and have been obtained by the same method as 
described above for the shear viscosity. The values of the several transport coefficients in 
the steady state are indicated. 

It is observed that the three Navier-Stokes transport coefficients are monotonically de¬ 
creasing functions of the speed parameter A* for the two values of the restitution coefficient 
a considered in the figures. A similar behavior was found for other values of a. The coeffi¬ 
cient of diffusive heat conductivity //, becomes even negative for large enough values of A*. 
Notice that this does not seem to imply the violation of any fundamental physical law or be 
incompatible with any physical symmetry. Nevertheless, it is quite possible that the exact 
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FIG. 3: (Color online) Adimensionalized shear viscosity rj of a two-dimensional system as a function 
of the dimensionless speed parameter A*. The (red) dashed curve is for a = 0.8 and the (black) 
solid line is for a = 0.9. The (blue) dots indicate the steady state values in each case. 



FIG. 4: (Color online) Adimensionalized (thermal) heat conductivity k of a two-dimensional system 
as a function of the dimensionless speed parameter A*. The (red) dashed curve is for a = 0.8 and 
the (black) solid line is for a = 0.9. The (blue) dots indicate the steady state values in each case. 


13 













0.4 



FIG. 5: (Color online) Adimensionalized diffusive heat conductivity JI as a function of the dimen¬ 
sionless speed parameter A* for a two dimensional system. The (red) dashed curve is for a = 0.8 
and the (black) solid line is for a = 0.9. The (blue) dots indicate the steady state values in each 
case. 

value of A* at which the change in sign of /.t occurs be a consequence of the approximations 
made and, in particular, of the first Sonine approximation. When the prediction for /i in 
this approximation is rather small, it is evident that higher order corrections might become 
relevant. 

IV. TRANSPORT COEFFICIENTS IN THE HOMOGENEOUS STEADY STATE 

A particularly relevant state of the confined granular gas is the homogeneous steady state. 
Stationarity of the temperature implies that the rate of change of the temperature vanishes, 

i.e. it is 

(‘"’(AiJ-O. (36) 
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Then, particularization of Eq. (IC3[) for the steady state yields 

da 2 


Ci, st — 2A 


St 


dA* 


Ci(A 


st) 


A*=A* 


X 


1 8x(a;,) 


+ Xi(a:,) 


4(^2 ,st + 1) — A* t 


<9a 2 
<9 A* 


a*=a* 




( 37 ) 


d(d + 2) 

Similarly, particularization of Eqs. (128]) - (130]) for the steady state lead to explicit expressions 
for the Navier-Stokes transport coefficients for this state, 

2 5 / 2 7T^ 


Vst = 


V T,\ A *> a 2 ,st), 


(d + 2)T (d/2) 71 


(38) 


^St 


2 5 / 2 (d-l)vr^ 
d(d + 2)T(d/2) 


1 + 2 a 2 st — 


A 


St 


da 2 


dA* 


a*=a: 


x 


02,st) + 


fJ'st 



1 -1 


A*=A* 


(39) 


(40) 


d(d + 2)T (d/2) ^(A* t , a 2 , s t) 

The frequencies 17^ and = 77^ are dehned in Eqs. (I34p . Moreover, the calculation of a 2 (A*) 
for A* in the vicinity of A* t can be carried out in an efficient and quite accurate way by 
noting that near A* st it is |<9a 2 /<9A*| -C 1 (see, for instance, Fig. 9 in Ref. 27]). Then, near 
the steady state, Eq. (jA13j) yields 

Bq + 4A 0 


a 2 


(41) 


B i + 4(A 0 + A\) 

The expressions of A 0 , Ai, B 0 , and B i are given in Eqs. flAl4l) - flAl7jl . 

It is now a simple task to evaluate the Euler and Navier-Stokes transport coefficients in 
the steady state. They are plotted in Figs. [6]|9] as a function of the coefficient of normal 
restitution a. The four coefficients present a clear dependence with the inelasticity. The 
Euler transport coefficient is a monotonic decreasing function of the coefficient of normal 
restitution, while the dimensionless shear viscosity monotonically increases with the value 
of a. This latter behavior is consistent with molecular dynamics simulation results reported 
both for dilute 29] and moderately dense systems 22]. Moreover, the dependence of the 
viscosity on the coefficient of normal restitution is clearly nonlinear, again in agreement 
with the simulations for dilute systems. It is worth to remind that in (non-conhned) dilute 
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FIG. 6: Dimensionless Euler transport coefficient of the two-dimensional confined granular gas in 
the homogeneous steady state Cl,st, as a function of the coefficient of normal restitution a. 


, 3 


the viscosity decreases as the coeffi- 


granular gases of smooth inelastic hard spheres 
cient of normal restitution increases, and that in a stochastic thermostat model it has been 

32]. On the other hand, the 


found to be a non-monotonic function of the inelasticity jill . 
dependence on the restitution coefficients of the two transport coefficients associated with 
the heat flux is not monotonic in the homogeneous steady state of the model discussed here, 
exhibiting both a minimum. Moreover, the coefficient of diffusive heat conductivity /z s t is 
negative in the whole range of values of a, while it is always positive in a dilute non-confined 
gas of inelastic hard spheres or disks. Notice that the dependence of the steady transport 
coefficients on the inelasticity of the system is quite strong. In particular, the (thermal) heat 
conductivity for a = 0.8 is about 25% smaller than its elastic limit value. 

In any case, when interpreting the results in Figs. [MHl it must be kept in mind that the 
bare transport coefficients have been scaled with a function of the temperature of the steady 
state, and that this temperature is in turn a function of the coefficient of normal restitution 
(and the velocity parameter of the model A). As a consequence, it is not possible to deduce 
the general expressions of the transport coefficients of the model to Navier-Stokes order from 
their form in the steady homogeneous state. 
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FIG. 7: Dimensionless shear viscosity of the two-dimensional confined granular gas in the homo¬ 
geneous steady state fj± st , as a function of the coefficient of normal restitution a. 



FIG. 8: Dimensionless (thermal) heat conductivity of the two-dimensional confined gas in the 
homogeneous steady state Tc s t, as a function of the coefficient of normal restitution a. 

V. DISCUSSION AND CONCLUSIONS 

The objective of this work has been to derive the hydrodynamic equations to Navier- 
Stokes order for a model of confined granular gas from an underlying kinetic theory, with all 
the parameters given explicitly. For clarification and context, the following comment must 
be taken into account. When applying the Chapman-Enskog procedure, the distribution 
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FIG. 9: Dimensionless diffusive heat conductivity of the two-dimensional confined granular gas in 
the homogeneous steady state Ji st , as a function of the coefficient of normal restitution a. 

function has been computed np to first order in the gradients of the hydrodynamic fields 
density, flow velocity, and granular temperature. Consequently, the heat and momentum 
fluxes are also determined to the same order. Since they occur as divergences in the balance 
equations, they lead to terms of second order in the gradients in those equations, what is 
usually referred to as the Navier-Stokes approximation for the fluxes. Also the rate of change 
of the temperature ( has been computed to first order, but it appears without any gradient 
operator in front of it in the balance equation for the energy. It follows that consistency 
of the Navier-Stokes order would require, in principle, computing ( up to second order in 
gradients, i.e. going an order further in the Chapman-Enskog expansion of the distribution 
function, the Burnett order. Such a calculation is rather involved and lengthy. The only case 
in which we are aware that second order contributions from ( have been analyzed deals with 
the linear contributions in a low density gas of smooth inelastic hard spheres Q. There, it 
is found than the terms are very small as compared with the similar ones arising from the 
fluxes. A similar behavior is likely to occur here with all the second order in the gradients 
contributions from (. 

The derived hydrodynamic equations are general, in the sense of having no restriction 
with regards to the values of the coefficient of normal restitution a or the velocity parameter 
of the model A. In particular, they hold in principle arbitrarily far from the homogeneous 
steady state, as long as the system be near an homogeneous hydrodynamic state. In the 
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steady state, both parameters a and A determine the temperature of the system, that is 
not an arbitrary parameter anymore. Of course, the general hydrodynamic equations can 
be particularized for the steady state, as it has been actually done here, but it must be 
emphasized that the general form of the hydrodynamic equations can not be inferred from 
the equations derived for the steady situation. 

The shear viscosity of the model in the steady state has been measured by using event 
driven molecular dynamics simulations. The transport coefficient was obtained from the 

29]. In a dense system (ncr 2 = 0.4) it was found 


decay rate of the transverse current 
that a linear fit in (1 — a) gives 


22 


rj ~ 1.0512^ [1 - 0.28(1 - a)} 

An expansion of Eq. (j38l) in powers on a to first order gives 

17. 

1 ~ 1 - 64 (1 - “)' 


(42) 


(43) 


As expected, the low density theory developed here is not able to predict the prefactor 
in Eq. (l42lh In dense systems, the collisional contributions to momentum transfer, and 
consequently to the shear viscosity play a fundamental role, and those effects are neglected 
at the level of the Boltzmann equation. However, if the lowest order inelasticity correction 
is considered, the results obtained in this paper are in good agreement with the simulation 
results for dense systems. A similar conclusion was reached in Ref. 29]. 

A relevant and recurrent question when deriving hydrodynamic equations from kinetic 
theory is to determine the context in which the equations apply. The small parameter in the 
Chapman-Enskog expansion is the ratio of the mean free path relative to the wavelength of 
the variation of the hydrodynamic fields. The mean free path is independent of the time 
for the homogeneous hydrodynamics. Consequently, it seems sensible to conclude that the 
conditions for the Navier-Stokes order hydrodynamics of the model are the same as for 
usual granular gases of inelastic hard spheres or disks and also for elastic collisions, i.e. for 
sufficiently large space and time scales as compared with the mean free path and the inverse 
collision frequency. 
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Appendix A: Chapman-Enskog solution 


To zeroth order in the gradients, the Chapman-Enskog expansion leads to 
dt 0 ) f { 0 \v)= [ dUiTotv, vi)/ ( 0 ) (u)/ ( 0 ) (t;i), 


while the balance equations to this order become 


df’n = 0, df ] u = 0, df ] T = -£ (0) T. 
The lowest order rate of change of the temperture is given by 

2 ;^[/ (0) ,/ (0) ]. 


C ( 0 ) M) = 


(Al) 


(A2) 


(A3) 


n(r, t)T(r, t)d 

It is worth to stress that the macroscopic fields are not expanded in the Chapman-Enskog 
method, so that the zeroth order distribution function already provides the exact actual 
macroscopic fields. Using Eqs. (lA2j) . Eq. (1A1I) can be transformed into 


-C < 0 ) AA- = / dv l T 0 (v,v 1 )f ( 0 , (v)f m (v l ). 


(A4) 


Since f(r,v,t ) and, therefore, /^°^(r,u,t) are normal, the latter must have the scaled form 

/ <0) (r, v, t ) = (c, A*), (A5) 


where 


vo m = 


2 T(r,t) 


m 


1/2 


(A6) 


is a characteristic local thermal velocity. The dimensionless function /®* only depends on 
the temperature through the scaled velocity, 


ce ^m 

v 0 {r,t ) 


and the dimensionless speed parameter 


A* 


A 


v 0 {r,t) 


(AT) 


(A8) 
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As a consequence, Eq. (1A4D is equivalent to 


C (°) 

~Y 


3 -(V/ (0) )+ 


dV 


dA* 


= / dniT 0 (n,ni)/ ( 0 ) (n)/ ( 0 ) (ni). 


(A9) 


This equation formally coincides with t 
of the system in the low density limit 


le one describing the homogeneous hydrodynamics 


26). ft is important to realize that the zeroth or¬ 


der approximation in the Chapman-Enskog method is not a local version, both in space and 
time, of the distribution function of the homogeneous steady state eventually reached by the 


system 


22 


25], but a local distribution generated from that describing the time-dependent 


homogeneous hydrodynamics. This is a relevant general issue when dealing with hydrody¬ 


namics around a non-equilibrium state 


37]. In the second Sonine approximation, the zeroth 


order distribution function is approximated by 


/ (0) (r, v, t) ~ 


n 


(V /2 i>o)' 


ie- c2 [l + a 2 (A*)S< 2 '(c 2 )], 


where 


g(2)( c 2^ = d(d + 2) _ d + 2 c 2 + ^ 
8 2 2 


(A10) 


(All) 


Neglecting nonlinear in a 2 terms, substitution of Eq. (IA10D into Eq. (IA3) yields 


C (0 \T) 


2 3/2 7r (d-l)/2 no .d-l t;o ( T ) 

r(d/2 )d 


a 




—) A* 2 
16 


(A12) 

An equation for a 2 (A*) is obtained by using this expression into the Boltzmann equation 
(IA9D . multiplication of the equation by c 4 , and later integration over the velocity c. If the 
quadratic terms in a 2 , as well as a term proportional to a 2 da 2 /d A*, are neglected, one gets 


26] 


da 2 

dA* 


4 4Ai + Bi 

-1- 

A* A 0 A* 


, 4 , B 0 

a 2 H-TT + 


A* A 0 A* ’ 


with the several coefficients given by 

Aq(«, A*) = (d + 2) 


1 — a 2 /7T\ 1 / 2 . . . 9 

-( - ckA - A 2 

2 \2J 


A 1 (a,A*) = 


(d + 2) 
16 


3(1 — a 2 ) 


+ A 


*2 


(A13) 

(A 14) 
(A 15) 


B 0 (a, A’) = (2ir) 1/2 (l + 2d + 3a 2 + 4A’ 2 ) aA - - 3 + 4A' 4 + a 2 + 2a 4 

—2d (1 - a 2 - 2A’ 2 ) + 2A* 2 (l + 6o 2 ) , (A16) 
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/ jr \ 1/2 I 

B^a, A*) = [2 - 2d(l - a) + 7a + 3a 3 ] A* - — {85 + 4A* 4 - 18(3 + 2a 2 )A 


- (32 + 87a + 30a 3 ) a - 2d [6A* 2 - (1 + a) (31 - 15a)] } 


(A17) 


The normal solution of Eq. (IA13D has been analyzed in Ref. 26], by solving numerically 
the differential equation for different initial conditions, and identifying the common part of 
all the generated solutions. Independently of the approximation used to calculate it, the 
solution of Eq. (1A9D is isotropic in velocity space, i.e. it only depends on V. It follows that 
the lowest order pressure tensor and heat flux are 


P l0) (r,«) =p(r,t)\, J„ (0) (r,i) = 0, 


(A 18) 


where I is the unit pressure tensor in d dimensions, and p(r, t ) is the hydrodynamic pressure 


P{r,t) = n{r,t)T(r,t). 

To first order in e, the expansion of the Boltzmann equation gives 

d/(°> 

dr 


(A19) 


a t (0) / (1> + £/ (1 > = /(»> 


V 


(A20) 


with the linear operator C defined by 


= - I dviT 0 (v,Vi) [/ (0) (r-.D.«)/ (1, (r.D 1 ,i) + f m (r,v,t)f m (r,v l ,t)} . 

(A21) 

The macroscopic balance equations to this order are 


d['^n + u ■ Vn + n'V ■ u = 0, 


d ( t ]> u + u ■ Vu + ( mn) 1 Vp = 0, 


2 T. 


<9 t (1) T + U-VT+ —V ■ u = -C (1) T, 
d 

with the first order in e cooling rate being a linear functional of 

4 


c (1) - 


[/ (0) ,/ (1) ] - 


(A22) 

(A23) 

(A24) 

(A25) 


n(r, t)T(r, t)d 

By using eqs. (lA22|b(lA24D . and that and are both normal distributions, Eq. (1A20|) 
is seen to be equivalent to 


< 9 f (0) / (1) +£/ (1) 


dfW 

dT 


C (1) T = A VlnT + B- VIn n + C : Vu, 


(A26) 
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where 


A{V\n,T) = 


Td f {0) . v 5 . (v/ (o:) ) + ^ A ^/ (0) 


2 dV 


dA* 


B(V\n, T) = -V/ (0) - 


T dfW 
m dV 


C( VKT) = A( V /(«))-I 




aw 


(A27) 

(A28) 

(A29) 


Because of the presence of the term involving A* on the right hand side of Eq. (1A29 


, the 


tensor C is not traceless, contrary to what happens in a system of elastic particles [34) and 
also in a system of inelastic hard spheres or disks jfij]. The solution of the linear equation 
(1A26|) must have the form 


/ {1) (V|n,T) = A VlnT+B-VlnT + C : Vu. 


(A30) 


Consider the first order cooling rate C*- 1 - 1 given by Eq. flA25[) . It is a scalar and, therefore, 
the only nonvanishing contribution to it has the form 


c (1) = ov- 


u 


(A31) 


with the Euler transport coefficient £1 given by 

4 


Ci = 


nTd 2 


o;[/ (0) ,TrC], 


where TrC denotes the trace of the tensor C. When Eqs. 


(A32) 


, (1A30D . and (1 A3 111 are 


substituted into Eq. (IA26I) . equations for A, B, and C are found by equating coefficients of 
the various gradients of the hydrodynamic fields, 


-C mT ^~ T ^A + CA = A, 

-C | 0 U-c ( 0 ) r|® + £B = B, 

- C + CC- Cl I = c. 
s dT d T s 


(A33) 

(A34) 

(A35) 


Next, let us analyze the contribution to the fluxes of first order in the gradients. The 
expression for pressure tensor contribution can be expressed as 


P (1) =m J dVVVfV>{V)=m J dV D{V)f {1) (V), 


(A36) 
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where 


D(V) = m 



(A37) 


since /(°' gives, by construction, the correct exact value of the hydrodynamic fields and, 
consequently, the contribituions to them from i = 1, 2, • • ■, must vanish. Then, taking 
into account the isotropy of the tensors and that D is traceless, it follows that 


P (1 ) = j dV D (V)C(V) : V« = —T] 
with (Vit) + being the transposed of Vu and 


Vu + (Vw) + —- V ■ u\ 

Lb 


(A38) 


W dvD(y):C(n 


(A39) 


Proceeding in a similar way, it is seen that the heat flux to first order in gradients reads 


jW = -kVT - (iVn. (A40) 

The coefficients in this expression are given by 



dV S(V) ■ A(V), 

(A41) 


K=~ d J dVS(V)-B(V), 

(A42) 

where 


(A43) 

Define the frequencies 

_ / dV D(V) : CC{V) 

" v ~ f dV D(V) : C(V) ’ 

(A44) 


_ / dV S(V) ■ £A(V) 

Uk ~ fdVS(V)-A(V) ’ 

(A45) 


_ / dV S(V) ■ jCB(V) 
iy,J ' ~ fdV S(V)-B(V) ' 

(A46) 


By using these definitions and Eqs. (1 A33 f) - ( 1 A35 [) . it is easy to obtain first order differential 
equations obeyed by the transport coefficients, 


(C (0) rA - r, = ~ d2+ 1 d _ 2 I dVD(V) : C(V) = —nT, (A47) 
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c^A -,. + M 

s dT dT 


K = 


^ J <*VS(V) ■ A(V) 

(d + 2)nT 


2 m 


1 + 2a 2 + T 


da 5 


^(0)rp_^ _ 




/i + 




K 


= ± j J d VS { V). B{ V ) = 


dT 
(d+2)T 2 


a 2 - 


(A48) 


(A49) 


dT ^ J n nd J ' v ’ ' v ’ ' 2m 

Upon deriving the last equalities in each of the three above equations, the second Sonine 
approximation for f(°\ Eq. (1A10[1 has been employed. 


Appendix B: Evaluation of the frequencies 


Approximated expression for the functions u v , u K , and have been obtained by using 
a Sonine expansion, truncated to lowest order. Since the collision operator commutes with 
the rotation operator, it follows from Eqs. f1A33[) - (1A35[) that A and B must be isotropic 
functions of the velocity times S(V). Moreover, only the traceless part of C is needed, and 
it must be the product of an isotropic function times D(V). Then, to lowest order in a 
Sonine expansion, 


A(V) oc f M (V)S(V), B(V) oc f M (V)S(V), C(V) TrC oc f M (V)D(V), 


where /m is the Maxwellian distribution, 


Im(V) = n ( 1 exp 


mV 2 


,2nTJ 1 V 2 T 

With this approximation, 

fdVD(V) : C[f M {V) D(V)] /dVD 0 -(V)L[/ M (V)D 0 -(V)] 

JdVf M (V) D(V):D(V) (d - l)(d + 2 )nT 2 

/ dV S(V) • C[f M (V) S(V)] 2 m / dV S(V) ■ C[f M (V) S(V)] 


Vr, = 




(Bl) 


(B2) 


(B3) 


(B4) 


f dV /m(V)S(V) ■ S(V) d(d + 2)nT 3 

The evaluation of the integrals in the above expressions is straightforward, and it is facilitated 
by using symbolic computer programs. Since similar calculations have been reported many 
times in the literature, we merely report here the final results, 

l MT) 


\f2i{ d 2 1 na d ln 


Vrj 


d(d + 2 )T(d/ 2 ) 


(1 + a) (2d + 3 — 3a) (l - g) 


+V2^(d - 2a) A* “ 2 ( 1 + 3 ^ J A 


3 CLr. 


*2 


(B5) 
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= 


7 r 2 n<j d 1 w 0 (T) f (1 + a) [512 + 352c? — 96a(d + 8)] 


V2d(d + 2)r (d/2) 

(1 + a) [5d + 4 — 3(4 — d)a] 


2 6 


2 6 


&2 


— (f) 1/2 [2(1 - d) + (rf + s) a ] A" 

A* 


32(d + 8) — 3(4 — d)a,2 ,* 2 


2 5 


For A* = 0, the above expressions reduce to those obtained in Refs. 


(B6) 


jp| and [6]. 


Appendix C: The Euler transport coefficient 


The transport coefficient A is given by Eq. (IA32D . An equation for it can be obtained 
from Eq. (1 A3 511 . 


- + £ T lC - T d J^ (l d = Tr C = 


(Cl) 


Since Tr C must have vanishing velocity moments up to second degree, its lowest order Sonine 
approximations reads 


TrC ~ 


b 2 e~ 


7r d/2 cr a!-l t ,rf+l 


SM(c 2 ). 


(C 2 ) 


To determine the dimensionless coefficient b 2 , we substitute the above expression into Eq. 
CD. and afterwards multiply it by v 4 and integrate over v. After some lengthy but trivial 
algebra, a differential equation is obtained, 


^(o) A * db 2 


dA* 


o7 (°) , 8 x( A *) . , ,> / . .x ,, . * da 2 


3C + 


d(d + 2) 


+ 4:d(i(a 2 + 1 ) — d(\A 


dA* 


bn = —2A : 


da 2 

dA* 


(C3) 


ffere, 


c = - = 

Ci b 2 


A _4+l 

47r 2 

d^~ 


dci / dc 2 e~ 


[l + a 2 SW(c?)] S<=>(cl) 


X 


A* 2 Ci 2 7T 1 / 2 aAc 2 2 (1 — tt 2 )c? 2 

Fppy + 2T (^) 4T pp) 


(C4) 


and 


x(A*) = 


J 


dvv l C Tr C 


dv f dvi f^(r, v, t) Tr C(vi)T 0 (v, vi)(v 4 + v 


(C5) 
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The integrals in the above expressions of (~i and x can be carried out getting 


Cl 


Ty-[96 + 9a 2 - 3a 2 (32 + 3a 2 ) + A* 2 (64 + 30a 2 )l , 

2 J fd 2 r(d/2) 1 J 


(C6) 


d-1 

X = 2 n T ( 2 rf/2) (y/2 {30« 4 (32 - a 2 ) - 5(544 + 7a 2 ) - 4A* 2 (32 + 15a 2 ) 

-64 (d - l)a(16 + a 2 ) - 2d(992 + 17a 2 ) + 3a 2 [928 + 43a 2 
+12A* 2 (32 + 3a 2 ) + 10d(32 - a 2 )] + 6A* 2 (288 - 45a 2 + 64d + 6 da 2 )} 
+512\/vrA* [2 + 7a + 3a 3 - 2d(l - a)]) . (C7) 


Now, the dimensionless coefficient 6 2 is obtained by solving numerically the set of equations 
(1A13I) and (IC3j) . and identifying the hydrodynamic part of the solution as discussed in the 
main text for the shear viscosity. Afterwards, the Euler transport coefficient (~i follows by 
using Eq. flC4jl . 
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